This section is intended to introduce various aspects of the art, which may be associated with exemplary embodiments of the present invention. This discussion is believed to assist in providing a framework to facilitate a better understanding of particular aspects of the present invention. Accordingly, it should be understood that this section should be read in this light, and not necessarily as admissions of prior art.
An important goal of seismic prospecting is to accurately image subsurface structures commonly referred to as reflectors. Seismic prospecting is facilitated by obtaining raw seismic data during performance of a seismic survey. During a seismic survey, seismic energy is generated at ground level by, for example, a controlled explosion (or other form of source, such as vibrators), and delivered to the earth. Seismic waves are reflected from underground structures and are received by a number of sensors/receivers referred to as geophones. The seismic data received by the geophones is processed in an effort to create an accurate mapping of the underground environment. The processed data is then examined with a goal of identifying geological formations that may contain hydrocarbons.
Seismic energy that is transmitted in a relatively vertical direction into the earth is the most likely to be reflected by reflectors. Such energy provides meaningful information about subsurface structures. However, the seismic energy may be undesirably diffused by anomalies in acoustic impedance that routinely occur in the subsurface environment. Diffusion of seismic energy during a seismic survey may cause subsurface features to be incorrectly represented in the resulting seismic data.
Acoustic impedance is a measure of the ease with which seismic energy travels through a particular portion of the subsurface environment. Those of ordinary skill in the art will appreciate that acoustic impedance may be defined as a product of density and seismic velocity. Acoustic impedance is typically referred to by the symbol Z.
Seismic attenuation is a common place in the Earth model which leads to amplitude loss and phase distortion in the seismic data. Seismic data that are affected by attenuation without proper treatment may lead to poor imaging of subsurface structures within and below the attenuating media. In seismic processing, attenuation compensation can be either carried out in a pre-processing stage or during a final imaging stage.
Seismic waves attenuate for a variety of reasons as they travel in a subsurface environment. A quality metric (sometimes referred to a quality factor) Q is typically used to represent attenuation characteristics of underground formations. In general, Q is inversely proportional to seismic signal attenuation and may range from a value of zero to infinity. More specifically, Q is a dimensionless quality factor that is a ratio of the peak energy of a wave to the dissipated energy. As waves travel, they lose energy with distance and time due to spherical divergence and absorption. Such energy loss must be accounted for when restoring seismic amplitudes to perform fluid and lithological interpretations, such as amplitude versus offset (AVO) analysis. Structures with a relatively high Q value tend to transmit seismic waves with little attenuation. Structures that tend to attenuate seismic energy to a greater degree have lower Q values.
Q values associated with subsurface structures are used to mathematically alter seismic data values to more accurately represent structures in the subsurface environment. This process may be referred to as “Q migration” by those of ordinary skill in the art. During Q migration, a seismic data value representing travel of seismic energy through a subsurface structure having a relatively low Q value may be decreased in amplitude and narrowed in spectrum to a greater degree than a data value representing travel of seismic energy through a subsurface structure having a relatively high Q value. Altering the amplitude and phase of data associated with low Q values takes into account the larger signal attenuation that occurs when seismic energy travels through structures having a relatively low Q value.
A crude attenuation compensation can be done as a pre-processing step by applying a global inverse Q-filter (usually 1D) to seismic data (Bickel and Natarajan, 1985; Hargreaves and Calvert, 1991; Wang, 2002), where corrections of both amplitude and phase, or either one of the two, are possible.
Compensating for attenuation during imaging allows for a more accurate and physical correction based on a Q model that can vary spatially. The Q model is usually provided by the users using Q tomography or some other field measurements. The Q compensation can be implemented in Kirchhoff depth migration (KDMIG) (Traynin et al., 2008), one-wave wave equation migration (WEM) (Mittet et al., 1995), and reverse time migration (RTM) (Deng and McMechan, 2008), all aiming for enhanced image quality by correcting reflector depth and improving illumination.
Q compensation in KDMIG and WEM is usually stable due to the fact that both methods have direct control over frequencies. On the contrary, time-domain RTM, which does not usually have control over frequencies, often suffers from exponential growth of energy in the high frequencies during Q compensation.
Conventionally, Q-compensated migration can be implemented in a few ways, which are discussed below.
Kirchhoff Depth Migration. At each imaging point, data are filtered (F) prior to applying the imaging condition using the following equation:
                              F          ⁡                      (                                          t                *                            ⁡                              (                t                )                                      )                          =                  exp          ⁢                      {                                                            ∓                  ω                                ⁢                                                                            t                      *                                        ⁡                                          (                      t                      )                                                        2                                            ±                              i                ⁢                                                                  ⁢                ω                ⁢                                                                            t                      *                                        ⁡                                          (                      t                      )                                                        π                                ⁢                ln                ⁢                                  ω                                      ω                    0                                                                        }                                              (        1        )            where t is the real-valued travel time, ω is angular frequency, ω0 is angular frequency, and t*(t) is the attenuated traveltime defined by:
                                          t            *                    ⁡                      (            t            )                          =                              ∫                                                                                Q                                          -                      1                                                        ⁡                                      (                    s                    )                                                                    v                  ⁡                                      (                    s                    )                                                              ⁢              ds                                =                      t                                          Q                eff                            ⁡                              (                t                )                                                                        (        2        )            which is based on integrating velocity v and Q along the rays s, wherein t is travel time and Qeff is the effective attenuation quality factor; but here, a ray-shooting strategy from the surface is used, as in standard Kirchhoff migration. (Cavalca et al., 2013).
Wave Equation Migration. Wave equation migration (Valenciano, et al. 2011) uses the following formula to downward continue wavefields:P(z+Δz;ω,kx,ky)=P(z;ω,kx,ky)exp(±ikzΔz)  (3)where P is pressure wavefield, ω is angular frequency, and the k's are the wavenumbers in the Cartesian x, y, and z directions. For visco-acoustic media, the normalized vertical wavenumber is given by:
                                          s            zQ                    =                                                                      (                                      1                    +                                          i                                              2                        ⁢                        Q                                                                              )                                2                            -                              s                x                2                            -                              s                y                2                                                    ⁢                                  ⁢        where        ⁢                                  ⁢                                            s              zQ                        =                                                            v                  ⁡                                      (                    ω                    )                                                  ω                            ⁢                              k                zQ                                              ,                                    s              zQx                        =                                                            v                  ⁡                                      (                    ω                    )                                                  ω                            ⁢                              k                x                                              ,                                    and              ⁢                                                          ⁢                              s                zQy                                      =                                                            v                  ⁡                                      (                    ω                    )                                                  ω                            ⁢                                                k                  y                                .                                                                        (        4        )            The frequency dependent velocity is
      v    ⁡          (      ω      )        =                              v          ⁡                      (                          ω              r                        )                          ⁡                  [                      1            -                                          1                                  π                  ⁢                                                                          ⁢                                      Q                    ⁡                                          (                                              ω                        r                                            )                                                                                  ⁢                              ln                ⁡                                  (                                      ω                                          ω                      r                                                        )                                                              ]                            -        1              .  Here Q(ωr) and v(ωr) are respectively the quality factor and phase velocity at a reference frequency.
RTM. Zhang (2010) has formulated a two-way wave equation for visco-acoustic media to be
                                          P            tt                    +                                    Φ              Q                        ⁢                          P              t                                      =                                            Φ              2                        ⁢            P                    =          0                                    (        5        )            where the subscripts t and tt denotes first- and second-order time derivatives, Φ is a pseudo-differential operator in the space domain defined by
      Φ    =                  (                              √                          (                                                -                                      v                    0                    2                                                  ⁢                                  ∇                  2                                            )                                            ω            0            γ                          )                    1                  1          -          γ                      ,v0 is the velocity defined at a reference frequency ω0, and ∇2 is the Laplacian operator and
  γ  =            1      π        ⁢                            tan                      -            1                          ⁡                  (                      1            Q                    )                    .      
Equation 5 can be modified to model amplitude-loss only as
                                          P            tt                    +                                    Φ              Q                        ⁢                          P              t                                -                                    v              2                        ⁢                                          ∇                2                            ⁢              P                                      =        0                            (        6        )            and phase distortion only asPtt+Φ2P=0  (7)The formulation in equation 7 was later extended to VTI (vertical transversely isotropic) media by Suh et al. (2012). Xie et al. (2015) proposed to implement attenuation compensating in TTI (tilted transversely isotropic) RTM based on the following visco-acoustic wave equation in TTI anisotropy:
                                                        (                              i                ⁢                                                                  ⁢                ω                            )                                      2              -                              2                ⁢                γ                                              ⁢                      (                                                            P                                                                              R                                                      )                          =                                                            v                0                2                            ⁢                              cos                2                            ⁢                                                π                  ⁢                                                                          ⁢                  γ                                2                                                    ω              0                              2                ⁢                γ                                              ⁢                      (                                                                                1                    +                                          2                      ⁢                                                                                          ⁢                      ɛ                                                                                                                                  1                      +                                              2                        ⁢                        δ                                                                                                                                                                                    1                      +                                              2                        ⁢                        δ                                                                                                              1                                                      )                    ⁢                      (                                                                                                      G                      xx                                        +                                          G                      yy                                                                                        0                                                                              0                                                                      G                    zz                                                                        )                    ⁢                      (                                                            P                                                                              R                                                      )                                              (        8        )            where P and R are the pressure and auxiliary wavefields, ε and δ are the Thomsen anisotropy parameters, Gxx, Gyy and Gzz are rotated differential operators to account for tilting of symmetry axis.